Time-resolved RNA signatures of CD4+ T cells in Parkinson’s disease

Parkinson’s disease (PD) emerges as a complex, multifactorial disease. While there is increasing evidence that dysregulated T cells play a central role in PD pathogenesis, elucidation of the pathomechanical changes in related signaling is still in its beginnings. We employed time-resolved RNA expression upon the activation of peripheral CD4+ T cells to track and functionally relate changes on cellular signaling in representative cases of patients at different stages of PD. While only few miRNAs showed time-course related expression changes in PD, we identified groups of genes with significantly altered expression for each different time window. Towards a further understanding of the functional consequences, we highlighted pathways with decreased or increased activity in PD, including the most prominent altered IL-17 pathway. Flow cytometric analyses showed not only an increased prevalence of Th17 cells but also a specific subtype of IL-17 producing γδ-T cells, indicating a previously unknown role in PD pathogenesis.


INTRODUCTION
Parkinson's disease (PD) is one of the most common neurodegenerative disorders with a still increasing worldwide incidence [1][2][3]. Cardinal PD motor symptoms include tremor, bradykinesia, rigidity and postural instability, which are connected to the loss of dopaminergic neurons of the substantia nigra pars compacta [4]. The pathological formation of α-synuclein aggregates and Lewy bodies has been associated with neuronal degeneration [4][5][6].
However, it is becoming increasingly evident that factors outside the brain are also significant contributors to PD pathogenesis. Deregulated immune cells could play a central role in the according scenarios [7][8][9]. Especially CD4+ T cells are of major importance for a coordinated immune function and appear to be deregulated in context of PD as they invade the brain and likely induce an autoimmune reaction [10][11][12]. Since functional T cell deregulation seems to emerge long time before the manifestation of PD cardinal motor symptoms, T cell changes also bear a great potential of being utilized as biomarkers for early diagnostics [12,13]. Additionally, T cell changes offer themselves as a starting point for the development of novel (immune) therapeutic strategies [8,14]. However, the pathomechanistic deregulation of T cell pathways awaits further clarification [8,15,16].
Time-course RNA expression data are particularly suitable to track and functionally relate changes on cellular signaling pathways [17,18]. Thus, to decipher deregulated T cell signaling in context with PD, we analyzed time-resolved RNA expression dynamics upon the activation of peripheral CD4+ T cells.
Based on the examination of five representative PD cases at different stages of disease, our study highlights comprehensive RNA expressional deregulations and identifies pathways related to changes on T cell functionalities in PD. Our data point to a mitochondrial dysfunction and a deregulation of DNA repair and are indicative for an increase of various autoimmune features and also provide first evidence for a yet unknown role of IL-17 producing gamma delta T cells in PD pathogenesis.

Experimental setup and quality control
To uncover deregulated T cell signaling in Parkinson's disease (PD), we analyzed peripheral CD4+ T cells from five PD patients in an age range between 53-86 and at different stages of disease i.e., stage 4, 4, 2.5, 1, and 1 according to the Hoehn & Yahr scale. Patient 5 was analyzed immediately upon diagnosis. As control we used five blood controls (healthy controls (HCs)) from healthy donors, aged 53-69. Since the initial 24 h time frame of the activation is decisive for a proper T cell immune function, including shifts of central signaling pathways and the transition from a resting to a proliferative cell stage [17,19,20], we induced the activation of the isolated CD4+ T cells by in vitro stimulation. Time-course sampling was done at different time-points (0, 2, 4, 8, 12, and 24 h), resulting in a total of 60 RNA samples (PD (n = 30); HC (n = 30)). These RNA samples were used for high-throughput analysis of time-resolved RNA expression profiles. An overview of the experimental setup is given in Fig. 1A. RNA integrity number (RIN) values of the extracted total RNA ranged from 7.1 to 9.1 indicating a high RNA quality for all time-course samples. Since CD69, IFIT3, and NME1 are reliable markers for the early 24 h T cell activation phase [21], we used them to verify the effective induction of the T cell activation process. Comparability between the two groups (PD patients and HCs) was confirmed by overlapping CD69, IFIT3, and NME1 mRNA expression patterns ( Fig. 1B-D).
Detection of differential transcriptomics during the time-course of T cell activation Analyses on the time-course transcriptional data by data clustering and multidimensional scaling identified distinct expressional changes between the PD and HC group. In detail, we calculated the pairwise Euclidean distance between the samples, across all genes and across all the time-points. We then applied the cmdscale function of the R base library to perform a classical multidimensional scaling with two dimensions (Parameters: eig = TRUE, k = 2). The combined time-course data allowed a clear separation between PD patients and HCs ( Supplementary Fig. 1A, horizontal axis). While patients at an early stage of PD (P4 and P5) mapped closer to HCs, patients with a progressive stage mapped at a greater distance to HCs, indicating of a possible link between clustering and severity or duration of PD. However, due to the small number of patients and controls, our study focuses on the separation between controls and PD patients considering both, early and late stages of the Parkinson´s disease.
A total of n = 30 158 genes were analyzed by our time-resolved RNA expression analyses. We found 535 genes with a median log 2 fold change of ≥0.5 or ≤−0.5 during the time-course. Out of the resulting 535 genes, 454 genes additionally showed a median log 2 fold change of ≥0.5 or ≤−0.5 for the comparison between PD patients and HCs. For the same comparison we found 175 deregulated genes with significantly different expression (p-values adjusted by Benjamini & Hochberg (adj. p-value ≤ 0.05) ( Fig. 2A).
We next addressed the question how many genes were differentially expressed in different time windows between PD patients and HCs. While expression changes were found throughout the entire time course, the significant changes were mainly attributed to the early and intermediary phases as detailed in the following. Likewise, multidimensional scaling under itemization of the different time-points, specified a clear distinction between PD patients and HCs, particularly for the early hours (0-4 h) of the T cell activation time-course ( Supplementary Fig. 1B).
Out of the 454 genes differentially expressed between PD and HCs with a median log 2 fold change of ≥0.5 or ≤ −0.5, there were 41 genes exclusively within the early time window between 0-2 h, 43 genes exclusively within the intermediary phase between 4-8 h, and 42 genes exclusively in the late phase between 12-24 h. There was an overlap of 102 genes that were shared between the early and intermediary phases, 38 genes shared between the intermediary and late phases and 167 genes that were found throughout all phases of the T cell activation course (Fig. 2B). Out of the 175 genes with a significantly altered expression, there were 72 genes specifically for the early time window and 32 genes specifically for the intermediary time window. There were no genes shared between the intermediary and late phases. Likewise, no genes were shared between the early and late phases. A single gene (HCAR3) was found in all phases of the T cell activation course (Fig. 2C). As exemplified for VCAN1 in the sections below, not all of the genes showed uniform directions of deregulations throughout the entire time-course.
Overall, we identified 172 core deregulated genes with both, a significantly altered expression (adj. p-value ≤ 0.05) and a median log 2 fold change of ≥0.5 or ≤ −0.5 for the comparison between PD samples and HCs ( Fig. 2A). In the following we address the different time points and provide examples of top deregulated genes (Fig. 3). Out of the core deregulated genes there were 117 genes with an increased RNA expression in PD at the 0 h timepoint. Some of the most extensively increased genes at this timepoint were IL12B, IL12A and TNFRSF8. TGFBI, ADAMTS10 and PACSIN1 were amongst the top genes with a decreased RNA expression in PD. At time points 2 h and 4 h, there were 53 and 117 of the core deregulated genes, respectively. Among the most extensively increased genes were IL12B, IL12A and CCL18 at the 2 h time-point and IL12A, CXCL9 and IL2RA at the 4 h time-point. Examples for the top genes with a decreased RNA expression were PACSIN1, PKMYT1 and RAD51 at the 2 h and PACSIN1, LOC100505585 and BRCA1 at the 4 h time-point, respectively. As for the later time points there were no deregulated genes with a significantly altered expression and a median log 2 fold change of ≥0.5 or ≤−0.5 at the 8 h and 12 h. However, at the 24 h time-point HCAR3 showed an increased expression meeting these criteria.

Specification of deregulated T cell gene expression signatures in context with PD
In the following, we focus on those core deregulated genes (i.e., genes with an adj. p-value ≤0.05 and median log 2 fold change of ≥0.5 or ≤−0.5 between PD and HC) that also show a significant deregulation at more than one time-point.
In total we identified 29 core deregulated genes with significantly increased RNA expression at least at two timepoints during the overall T cell activation time-course (Fig. 4A). Corresponding genes encode among others the chemokine ligand CCL18, the growth factor and type 2 cytokine AREG, the mitochondrial monoamine oxidase MAOA, the hydroxycarboxylic acid receptor HCAR3 and the TNF receptor TNFRSF8 (  A Following T cell activation a total of 30,158 genes was analyzed for time-course RNA expression changes. There were 535 genes with a median log 2 fold change (FC) ≥0.5 or ≤−0.5 during the time-course. The comparison between the PD and HC groups showed 454 genes with a median log 2 FC ≥0.5 or ≤−0.5 and 175 genes with a significantly altered expression (p-values ≤ 0.05 adjusted by Benjamini & Hochberg). The 172 genes that showed both a significantly altered expression and a median log 2 FC ≥0.5 or ≤−0.5 are referred to as core deregulated genes. B, C The Venn diagrams provide a numerical summary of genes that were deregulated at different time windows upon T cell activation. The diagram (B) includes the genes with a median log 2 FC ≥0.5 or ≤−0.5 and the diagram (C) the genes with an adjusted p-values ≤ 0.05.
There were 46 of the core deregulated genes with a decreased RNA expression at least at two time-points in PD samples during the T cell activation time-course (Fig. 5A). These encode among others the extracellular protease ADAMTS10, the nuclear phosphoprotein BRCA1, the uncharacterized gene LOC100505585 and the replication initiation factor MCM10 ( A single mRNA (extracellular matrix regulator VCAN1) revealed an alternation of increased RNA expression at 0 h (with a log 2 fold change of 1.20 and an adj. p-value of 3.13 × 10 −2 ) and a decreased RNA expression at 4 h (with a log 2 fold change of −0.66 and an adj. p-value of 4.86 × 10 −2 ) in PD samples (Fig. 5F, G).
Notably and as summarized in Supplementary table 2, the above-mentioned deregulated genes can be related to PD, neurodegeneration, and autoimmune diseases and are currently being under investigation for therapeutic use.
Specification of deregulated miRNA expression in context with PD We next performed time-resolved expression analyses of the miRNome. In contrast to the transcriptome analyses we identified only very few miRNAs, according to the above applied criteria (i.e., time-course expression changes and deregulation in PD patients with significant differences at more than one time-point) (Fig. 6).
Hsa-miR-132-3p showed significantly increased expression at 8 h and at 24 h in PD samples. The according median log 2 fold changes between PD and HC samples were 0.59 and 0.58, corresponding FDR adjusted p-values were 1.66 × 10 −2 and 3.26 × 10 −2 , respectively. Hsa-miR-223-3p showed significantly increased expression at 0 h and between 8 h-24 h in PD samples. MiRNAs usually exert rather moderate effects as frequently shown by decreased mRNA levels of their respective target genes [22,23]. Putative targets of the above-described miRNAs were identified by inverse correlations of the according miRNA and mRNA time-course changes. For each time-point, median log 2 FCs were determined for the comparison between PD and HC groups. An inverse Pearson´s correlation coefficient (PCC) ≤−0.5 was chosen for the matching of resulting miRNA and mRNA timecourse changes. As an additional criterion, a median log 2 FC decrease level of at least −0.3 was chosen for putative targets of the miRNAs that showed an increased expression in context with PD. A median log 2 FC of ≥0.3 was chosen for putative targets of miRNAs with a decreased expression in PD. Corresponding analyses for hsa-miR-132-3p and hsa-miR-223-3p identified 1,525 and 862 target genes, respectively. Analyses for hsa-miR-4730, hsa-miR-155-5p and hsa-miR-762 identified 1,293, 1,747 and 1,101 genes, respectively, which are subject to direct or secondary miRNA regulatory effects. To identify direct miRNA targets we included information on strong experimental target validation as determined by miRTargetLink 2.0 [24] (Supplementary table 3). Amongst others, we identified CRK as target of miR-132-3p [25] with a link to T cell migration [26], and TP53INP1 as target of miR-155-5p [27][28][29][30][31][32][33] with a link to cellular stress response [34].

Increase of IL-17 signaling, Th17 and IL-17 producing gamma delta T cells in PD
Numerous of the genes that we found with an increased RNA expression pattern in connection with PD, could be assigned to the KEGG category "IL-17 signaling pathway" (as already mentioned above). Representative genes such as CCL20, IL12A, IL12B, IL17D, IL17F, IL21, IL22, IL6, PTGS2 and TNFRSR8 showed distinct mRNA deregulation during the time-course (Fig. 8A). RORC and BATF are considered as crucial transcription factors for the induction of IL-17 signaling and Th17 development [39,40]. Examination of the corresponding time-course expression data indicated increased mRNA levels of these key transcription factors as part of the T cell activation process (Fig. 8B, C).

DISCUSSION
Despite first evidence for an involvement of peripheral T lymphocytes in PD pathogenesis and their potential utility for novel therapeutic strategies, there is still a lack of studies that specify relevant deregulation of T cell pathways in context with PD [8,15]. To decipher deregulated T cell signaling in PD, we analyzed time-resolved RNA expression patterns (0-24 h) following the in vitro activation of peripheral CD4+ T cells.
When comparing the T cell activation time-course data of our PD patients to HCs, we identified an up-regulation for miR-223-3p. An increased prevalence of this miRNA has also been detected within the blood sera of PD patients [41] and a connection to brain autoimmune diseases has several times been reported [42,43].
Additionally, we found miR-155-5p with a significant downregulation in the CD4+ T cell time-course data of PD patients. This miRNA has, however, formerly been reported to have an increased prevalence in the blood cells of PD mouse models [44]. Our findings of a significantly reduced miR-155-5p expression were mainly attributable to the PD patients 1-3. It is legitimate to speculate that the detected effect may be related to treatment with levodopa [45], which was given to the patients 1-3. Admittingly, it is highly challenging to pinpoint drug related effects in studies with a small number of patients that received different medical treatments. We do, however, feel that the overall RNA data in our study does not reflect the treatment schemas: While we found RNA deregulations common to all patients, each patient was treated by different drugs and drug combinations.
We also identified significant deregulation of miR-762. This miRNA has been described as a regulator of energy metabolic functions and the production of reactive oxygen species (ROS) in mitochondria [46]. Elevated production of ROS is a well-known reason for the emergence of DNA damages and is commonly associated with pathological changes in cell physiology [47]. Our finding of a miR-762 deregulation may be indicative for the according scenarios in CD4+ T cells of PD patients. We also observed deregulations of specific mRNAs, including BRCA1, MAOA, MCM10, or RAD51, which encode mitochondrial proteins and enzymes involved in maintaining the genomic integrity [48][49][50]. Additionally, our enrichment analysis of the transcripts that were decreased in PD identified functional links to synthesis-dependent strand annealing (SDSA) that is part of the DNA double-strand break repair pathway [51]. Our findings indicate that mitochondrial dysfunctions and a malfunction in DNA repair mechanisms likely contribute to T cell related PD pathogenesis. This assumption was further supported by analyses of a specific sub-group of alpha-synuclein specific memory T cells in PD [52]. Common genomic mutations, contributing to the development of PD, also relate to changes of mitochondrial functions and DNA repair mechanisms [53][54][55]. Deregulations of these processes are considered as central hallmarks of PD in neuronal cells [56]. However, the functional effects of genomic mutations may possibly not be restricted to neurons but may also have a bearing for pathological changes in other cell types, including the CD4+ T cells.
Further RNA expressional deregulation can be assigned to central immune and cytokine pathways including the IL-10, Jak-STAT and G alpha (i) signaling pathways. Shifts on corresponding pathways are commonly associated with pathological T cell functions and immune disorders [57][58][59]. As already mentioned in the results section, the deregulation of representative transcripts such as AREG, CCL18, HCAR3, and TNFRSF8, further supports the contribution of T cell autoimmune features to PD pathogenesis.
The IL-17 signaling pathway in particular appears to play a prominent role for the development of brain related autoimmune diseases [60]. Notably, animal models also provided evidence for a link between circulating IL-17 and cognitive impairments [61]. In our analyses we found strong evidence for the deregulation of IL-17 signaling. In detail, the subtype composition of corresponding peripheral CD4+ T cell samples revealed a higher incidence of cells that were assigned to the Th17 subtype. In line with our results, Sommer et al. showed an increased prevalence of IL-17-producing T cells in patients at early stages of PD [62]. Interleukin-1 (IL1) that is secreted by activated microglia in context with PD pathogenesis [63] is a critical cytokine for Th17 cell differentiation [64,65]. As shown for other neurodegenerative diseases [66] a crosstalk between the brain and periphery may provide a potential explanation for our findings of increased TH17 counts in PD peripheral blood samples. In addition, we found an increased prevalence of a specific subtype of IL-17 producing gamma delta T cells in PD patients. To our knowledge, the involvement of this cell type in PD pathogenesis has not been reported yet. IL-17 producing gamma delta T cells can be activated, independently of any T cell receptor stimulation, through cytokines such as IL1 [67] and are considered to play an important role for autoimmune diseases by amplifying Th17 responses [68][69][70]. It is conceivable that this cell type contributes to the development of an autoimmune reaction as part of the PD pathogenesis.
In summary, our time-resolved RNA expression analyses of in vitro activated CD4+ T cells uncovered comprehensive RNA Over-representation analysis (ORA) and protein interaction network analysis of core deregulated genes in PD. A, C Corresponding FDR adjusted p-values are shown for the ORA results of "KEGG" and "Reactome" databases. Protein interaction networks generated for genes with increased (B) and decreased (D) expression were exported from STRING database (v 11.5).
expressional deregulation in PD. Associated cellular processes point to a mitochondrial dysfunction and a disruption of DNA repair in CD4+ T cells of PD patients. Various connections to autoimmunity and to central PD features are highlighted. We found strong indications for the pathomechanistic relevance of IL-17 signaling and provide first evidence for the involvement of IL-17 producing γδ-T cells in PD. Our analyses refer to a relatively small number of exemplary cases. Extended analyses with larger patient groups will help to assess their potential for future clinical implications.

MATERIALS AND METHODS PD patients and healthy controls
Epidemiological studies show an increased prevalence of PD in connection with progressing age and with the male sex [1,71,72]. Based on this, the PD patients (n = 5) of this pilot study constitute a cohort of elderly male (age 53-85 yrs.; non-smokers), representative for different stages of disease development by the Hoehn and Yahr Scale [73,74]. Corresponding HCs (n = 4) were matched for age (age 53-63 yrs.; non-smokers) and gender. One HC was carried out twice (C1 and C1_2), representing independent replicates from the same healthy donor but from two different days of blood collection ranges by single data points. P-values for the comparison between PD and HC groups are indicated for the individual timepoints by asterisks (ns: not significant; *p-value ≤ 0.05; **p-value ≤ 0.01). D-G Flow cytometric analysis was done for CD4+ T cells, 6 h after CD4+ T cell activation. D Purity of the isolated CD4+ T cells was confirmed for HC and PD groups. E, F Significant increase of Th17 (IL-17+) and gamma delta (γδTCR+) subtypes were detected amongst the CD4+ T cells of the PD group. G Significant increase of IL-17 producing T cells was detected amongst the γδTCR + CD4+ T cells subpopulation (γδTCR + IL-17+) of the PD group. Median results of HC and PD groups are indicated by horizontal lines and total percentage ranges by individual data points. P-values for the comparison between PD and HC groups are specified.
(at a time distance of seven months). Basic subject data, including age, stage of disease, reported PD onset, clinical features, comorbidities and medications, are summarized in Supplementary table 1. Peripheral blood samples for subsequent T cell isolation (27 ml per donor) were drawn to lithium heparin containing collection tubes (S-Monovette, Sarstedt AG & Co. KG, Numbrecht, Germany). Written informed consents were obtained from all subjects (ethics approval ID: Ethical vote No. 213/08) and a previous immune reaction was excluded to the time-point of blood collection by the analysis of total blood counts from an additional aliquot (2 ml). CD4+ T cell isolation and collection of activated time-course samples CD4+ T cells were isolated from peripheral blood samples by density gradient centrifugation, followed by negative selection with the Human CD4+ T cell Isolation Kit (Miltenyi Biotech, Bergisch Gladbach, Germany). The isolated cells were suspended with RPMI 1640 medium (Life Technologies GmbH, Darmstadt, Germany; with 10% v/v heat inactivated fetal bovine serum (Biochrom GmbH, Berlin, Germany) and 1% v/v penicillin-streptomycin (100 U/ml)) and incubated in 25 mm flasks overnight. At the following day, the cells were seeded in a 96 well format (350,000 cells/well) and were in vitro activated by the human T cell activation/expansion kit (αCD2/αCD3/αCD28 MACSiBead particles, Miltenyi Biotec GmbH, Bergisch Gladbach, Germany). Cellular samples were collected from different wells at 0, 2, 4, 8, 12 and 24 h after activation for subsequent RNA extraction as formerly described [17]. Time-course collection resulted in a total of n = 30 samples from PD patients (n = 5) and controls (n = 4 + n = 1 (C1 replicate)), respectively.

RNA extraction and determination of time-resolved expression profiles
Cellular RNA was extracted, quality controlled, quantified and analyzed by microarray systems from Agilent Technologies (One-Color, Human SurePrint G3) as detailed in a previous publication [17]. For the timeresolved transcriptome analysis one sample (patient P5) yielded no data at the 4 h time-point.

Data processing
For both miRNA and transcriptome microarrays, raw expression values were processed using the limma R-package [75]. A background correction (method = normexp, offset = 16), quantile normalization, and log2 transformation were conducted. For the transcriptome data, a batch correction ('removeBatchEffect' method of the limma R-package) was applied to account for potential variations, due to different microarray batches.

Statistical analysis
For the statistical analysis of the expression values, we used a two-step approach. In the first step, we identified genes or miRNAs that showed a change in expression during the analyzed 24 h period after T cell activation. For these genes, we then compared expression differences between PD patents and HCs in the second step.
To determine the expression change of a given miRNA or gene in the analyzed time frame, we calculated their maximum deviation from the initial time point using the log2 fold-change for each sample individually. We then aggregated the scores for the PD and the HC groups, respectively, using the median. Finally, we select all genes or miRNAs that appeared with an absolute fold-change of at least 1.5 in one on the two analyzed groups.
To compare gene expression differences between PD patients and HCs, we calculated log 2 fold-changes between the median values in both groups and conducted a shrinkage t-test [76]. For all analyses, the p-values were FDR adjusted [77].

Enrichment analysis
All enrichment analyses in this study were performed using version 3.2 of the GeneTrail web service [35]. For the analysis of deregulated biological processes we applied over-representation analyses for each test set using all protein coding genes as a reference set. All resulting p-values were FDR adjusted [78].

Protein network analyses
Protein interaction networks were generated using the STRING database (v 11.5) [38], excluding text-mining from the selection of active interaction sources and hiding disconnected nodes within the networks.

DATA AVAILABILITY
The datasets of the current study are available on the GEO database repository. The super-series number of the study is #GSE202667. Accession numbers for timeresolved transcriptome profiles and time-resolved miRNA profiles are #GSE202665 and #GSE202666, respectively.